function [node]=vitUpdateLocalMapLocalization(Network,epochs,r)
%localMapLocalization computes the local 1 or 2 hop map for all its nodes
%coordinates for each node in the 'Network'. It intends to be faster than  
% the original localMapLocalization function. (Try to
%use this to replace 'neighborMap' and 'ccaLocation_hop')
% Network - input deployed network
% epochs - required number of training cycles
% r -radius
% node_k - integer of the node number
% T - elapsed time for CCA computation
% C_delta - difference of C_tranform from real coordinates. C_tranform is the 
% the resulting node coordinates in 2D
%
% D_dist_mean - the sum of the mean error between the original distance matrix of the input
% 'Network' and the distance matrix computed from the resulted node coordinates
%
% D_coordinates_mean - the mean error between the coordinates of the input
% 'Network' and the resulted node coordinates
%
% D_coordinates_median - the median error between the coordinates of the input
% 'Network' and the resulted node coordinates

THRESHOLD=0.5; %the max distance error average to be tolerated or we'd fail
% THRESHOLD=5; %used when measurement errors are big
r  %just to display what radius we are computing on the screen. 
N=size(Network,1);
D=sqrt(disteusq(Network,Network,'x'));
%added block on May 8th for getting a better noise level
D_measure=D;
tmp = 0.05*(randn(N)); %added 5% error of normal distribution
for ii=1:N
    for jj=ii+1:N
        tmp_sym(ii,jj)=tmp(ii,jj);
        tmp_sym(jj,ii)=tmp_sym(ii,jj);
    end
    tmp_sym(ii,ii)=0;
end
% ldist=ldist.*(1+tmp);
D_measure=D_measure.*(1+tmp_sym);
%added block on May 8th

cl_count=0;
for i=1:N %compute all the local maps
%get node_i's neighbor nodes within r;
[node(i).neighbors]=find_neighbors(D,r,i,1);
cl_count=cl_count+size(node(i).neighbors,2)-1; %trying to get network connectivity level
end
cl_count=cl_count/N; %network connectivity level
%D=sqrt(disteusq(Network,Network,'x')); 
%N=size(Network,1);
%%%%%%%%%%%%%%%%% this was used in rectangle 4x25 network %%%%%%%%%%%%
% if cl_count<=5.2 
%     t_level=15;
% end
% if (cl_count>5.2)&(cl_count<12.02) 
%     t_level=12 %select connectivity level of 12 as the cut for having h=1
% end
% if cl_count>=12.02
%     t_level=15 % increase the matrix size if network connectivity level is high
% end % manipulate the size of the local map. Experimental only. May be different for each configuration
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% if cl_count<12.02
%     t_level=12;
% else t_level=15;
% end % used for random 200 square and cshape_79nodes_grid network
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% t_level=20 %use this for rectangle grid networks also quite O.K. 
%%%%%%%%%%%%%%%%%%%%%%%%%used for cshape_79nodes_grid network
% t_level=12; %used for square grid. 
t_level=20; %used for random loop network; for cshape random

for i=1:N
    for j=1:N
        if D(i,j)<r
            D_hopDist(i,j)=D(i,j);
        else
            D_hopDist(i,j)=2*N*r;
        end
    end
end %get prepared to compute the shortest distance matrix for D

for k=1:N
    D_hopDist = min(D_hopDist,repmat(D_hopDist(:,k),[1 N])+repmat(D_hopDist(k,:),[N 1])); 
end %compute the shortest distance matrix using Floyed algorithm

% for k=1:63
%     for i = 1:63
%         for j = 1:63
%             D_hopDist2(i,j) = min(D_hopDist2(i,j), D_hopDist2(i,k)+D_hopDist2(k,j));
%         end
%     end
% end %compute the shortest distance matrix using Dijkas. Performance is bad in matlab

for i=1:N
    for j=1:N
        if D(i,j)<r
            D_hopCount(i,j)=1;
        else
            D_hopCount(i,j)=2*N;
        end
        if(i==j) D_hopCount(i,j)=0;
        end
    end
end %get prepared to compute the shortest hop matrix for D

for k=1:N
    D_hopCount = min(D_hopCount,repmat(D_hopCount(:,k),[1 N])+repmat(D_hopCount(k,:),[N 1])); 
end %compute the shortest hop matrix using Floyed algorithm

for node_k=1:N
    node_k
    tStart = cputime;
    while(size(node(node_k).neighbors,2)==1)%node_k has no connectivity
        node_k=node_k+1;
    end
    %tic;
    h=2;
    
    node(node_k).t_level=t_level;
    if (size(node(node_k).neighbors,2)>t_level) 
        h=1;
    end
    %added D_measure in the line below on May 8th 
    [local_d,node_index]=localDist(D_hopCount,D_hopDist,D_measure,r,node_k,h); %compute node_k's h hop distance matrix
    local_size=size(node_index,2);

    %node_index_merge=node_index;
    %local_d_merge=local_d;
    %local_size_merge = size(node_index_merge,2);
    %get the merged network real value from the deployed network
    local_network=Network(node_index,:); %use this line to replace the following three lines
%     local_network=Network(node_index(1),:);
%     for i=2:local_size
%         local_network=[local_network;Network(node_index(i),:)];
%     end %this loop is not needed. Simplify it.
    %compute the real distance from the deployed network node values
    local_distance_deployed=sqrt(disteusq(local_network,local_network,'x'));
    delta_distance=local_d-local_distance_deployed;
    node(node_k).Dmatrix_error_orig_mean = mean((mean(abs(delta_distance)))');

    node(node_k).local_d_merge=local_d;
    node(node_k).local_distance_deployed=local_distance_deployed;
    node(node_k).local_network=local_network;
    node(node_k).neighbors_merge=node_index;
    node(node_k).node_id=node_k;
    node(node_k).radius=r;

    %calculate the proximity of the matrix using SVD - not useful
    %[u,s,v]=svd(local_d_merge);
    %local_d_merge=u(:,1)*s(1,1)*v(:,1)'+u(:,2)*s(2,2)*v(:,2)'+u(:,3)*s(3,3)*v(:,3)'+u(:,4)*s(4,4)*v(:,4)';

    %calculate the local network value using cca
    pass=0;
    node(node_k).D_error_mean_compute=THRESHOLD; %set the threshold
    epochs_adjust=epochs;
    epochs_max=250;

    for jj=1:5
        %tic;
        local_network_c=cca(local_d,2,epochs_adjust,local_d);%computed local_network of node_k
        %T=toc;
        D_C = sqrt(disteusq(local_network_c,local_network_c,'x'));
       % D_dist_mean = mean((mean(abs(D_C-local_distance_deployed)))');
        D_dist_mean = mean((mean(abs(D_C-local_d)))');
        D_dist_mean=D_dist_mean/r;
        if(D_dist_mean<node(node_k).D_error_mean_compute)
           node(node_k).epochs=epochs_adjust;
           node(node_k).local_network_c=local_network_c;
           %node(node_k).local_map_compuTime=T;
           node(node_k).D_error_mean_compute=D_dist_mean;
        end
        if (D_dist_mean<0.03)&(local_size>8) %used 0.04 for low connectivity and 0.02 for high connectivity
            pass=1;
            %T=toc;
            T=cputime-tStart;
            break;
        end %if
    end% for j

    if (pass==0)&(D_dist_mean>0.05) %try some other epochs values 
        if (local_size>=10) 
            epochs_max=150;
        end
        if (local_size>=8)&(local_size<10)
            epochs_max=200;
        end 
            while (epochs_adjust<epochs_max)
                epochs_adjust=epochs_adjust+50; %increase the epochs
                if (local_size>19) 
                    mini_cycle=3; %can use different mini_cycle of big and small matrixes
                else mini_cycle=3; %used 5, 3 and other numbers, all work O.K.
                end
                
                for j=1:mini_cycle %give it another attempt
%                    tic;
                   local_network_c=cca(local_d,2,epochs_adjust,local_d);%computed local_network of node_k
                   %T=toc;
                   D_C = sqrt(disteusq(local_network_c,local_network_c,'x'));
%                  D_dist_mean = mean((mean(abs(D_C-local_distance_deployed)))');
                   D_dist_mean = mean((mean(abs(D_C-local_d)))');
                   D_dist_mean=D_dist_mean/r;
%                    
                   if((D_dist_mean<node(node_k).D_error_mean_compute))
                       node(node_k).epochs=epochs_adjust;
                       node(node_k).local_network_c=local_network_c;
                       %node(node_k).local_map_compuTime=T;
                       node(node_k).D_error_mean_compute=D_dist_mean;
                   end
                   
                   if (D_dist_mean<0.03)&(local_size>8)
                       pass=1;
                       %T=toc;
                       T=cputime-tStart;
                       break;
                   end
                end %for j
                if (pass==1)
                    break; %break out of 'while'
                end
            end %while
    end %if (pass==0)
        
    
    if ((pass==0)&(node(node_k).D_error_mean_compute<THRESHOLD))
%         T=toc;
        T=cputime-tStart;
        node(node_k).local_map_compuTime=T;
        T=0;
    end %we run all the epochs values. 

    if (pass==1)%D_dist_mean<0.04
%         node(node_k).epochs=epochs_adjust;
%         node(node_k).local_network_c=local_network_c;
        node(node_k).local_map_compuTime=T;
        T=0;
    end

    if ((pass==0)&(node(node_k).D_error_mean_compute==THRESHOLD))
        fprintf(2,'node %d failed to compute the cca\n', node_k);
    end %we can't get it for this node
    
    %debug
%     if((D_dist_mean>=0.02)&(epochs_adjust<250))
%         fprintf(2,'early quit of cca sycles\n');
%     end

    %check the accuracy of our results
    N_local=size(local_network,1);
    N1=1;
    N2=int16(N_local/3);
    N3=int16(N_local*2/3+1);
    Y(1,:)=node(node_k).local_network_c(N1,:);
    Y(2,:)=node(node_k).local_network_c(N2,:);
    Y(3,:)=node(node_k).local_network_c(N3,:);
    X(1,:)=local_network(N1,:);
    X(2,:)=local_network(N2,:);
    X(3,:)=local_network(N3,:);
    [d, Z, transform] = procrustes(X, Y);
    Cx=transform.c;
    for i=1:(N_local-size(transform.c,1))
          Cx=[Cx;transform.c(1,:)];
    end
    C_transform=transform.b*node(node_k).local_network_c * transform.T+Cx;

    % [Z,Cx]= mapTrans(X1,X2,X3,Y1,Y2,Y3,N);
    % C_transform=(Z*local_network_c'+Cx)';
    node(node_k).local_network_transform=C_transform;
    D_C = sqrt(disteusq(C_transform,C_transform,'x'));
    D_dist_mean = mean((mean(abs(D_C-local_distance_deployed)))');
    D_dist_mean=D_dist_mean/r;
    node(node_k).D_dist_mean_true=D_dist_mean;
    D_coordinates_mean=mean(abs(C_transform-local_network));
    D_coordinates_mean=D_coordinates_mean/r;
    D_coordinates_median=median(abs(C_transform-local_network));
    D_coordinates_median=D_coordinates_median/r;
    C_delta=C_transform - local_network;
    node(node_k).local_coordinates_error_mean=D_coordinates_mean;
    node(node_k).local_coordinates_error_median=D_coordinates_median;
end %for node_k

%%%%%%%%%%%%%%%%subfunctions%%%%%%%%%%%%%%%%%%%%%%%%%%
function [ldist,node_index]=localDist(D_hopCount,D_hopDist,D_measure,r,node_i,hop)
%added D_measure on May 8th
%localDist(D,r) takes a distance matrix D and range r to generate for  
%node 'i' 0<i<size(D) a local distance matrix ldist that includes the
%nodes j such that D(i,j)<=hop*r. All the unknown entries in ldist would be
%marked as 'NaN'. "node_index" is an array that holds the node index
%'j' of the node j in D that are selected and included in the ldist, in 
% the ascending order. Nodes in ldist also have their index in the 
% ascending order of the original index in D. 

% N=size(D,1);
% for i=1:N
%     for j=1:N
%         if D(i,j)<r
%             D_hopDist(i,j)=D(i,j);
%         else
%             D_hopDist(i,j)=2*N*r;
%         end
%     end
% end %get prepared to compute the shortest distance matrix for D
% 
% for k=1:N
%     D_hopDist = min(D_hopDist,repmat(D_hopDist(:,k),[1 N])+repmat(D_hopDist(k,:),[N 1])); 
% end %compute the shortest distance matrix using Floyed algorithm
% 
% % for k=1:63
% %     for i = 1:63
% %         for j = 1:63
% %             D_hopDist2(i,j) = min(D_hopDist2(i,j), D_hopDist2(i,k)+D_hopDist2(k,j));
% %         end
% %     end
% % end %compute the shortest distance matrix using Dijkas. Performance is bad in matlab
% 
% for i=1:N
%     for j=1:N
%         if D(i,j)<r
%             D_hopCount(i,j)=1;
%         else
%             D_hopCount(i,j)=2*N;
%         end
%         if(i==j) D_hopCount(i,j)=0;
%         end
%     end
% end %get prepared to compute the shortest hop matrix for D
% 
% for k=1:N
%     D_hopCount = min(D_hopCount,repmat(D_hopCount(:,k),[1 N])+repmat(D_hopCount(k,:),[N 1])); 
% end %compute the shortest hop matrix using Floyed algorithm
%all the above matrix computation should be moved out of this function into
%the main program so that it would not be repeated
% tic;
% tStart = cputime;
N=size(D_hopDist,1);
node_count=0;
for j=1:N
    if D_hopCount(node_i,j)<=hop %look for all the nodes that is within hop 
        node_count=node_count+1; %count the number of selected nodes
        node_index(node_count)=j; %save the original index of the selected node
%         for m=1:N
%             if ((m~=i)&(D(j,m)<=r))
%                 node_count=node_count+1;
%                 node_index(node_count)=m;
%             end
%         end
    end 
end %get all the nodes within the hop range.  
%node_index=sort(node_index);
ldist=zeros(node_count,node_count);
for i=1:node_count
    for j=i+1:node_count
        ldist(i,j)=D_measure(node_index(i),node_index(j)); %used D_measure on May 8th
        ldist(j,i)=ldist(i,j);
%         if (ldist(i,j)>r) 
%             ldist(i,j)=-1;
%             %ldist(i,j)=ceil(ldist(i,j)/r)*r;
%             %if (k==1) ldist(i,j)=(2*r+r)/2; end
%         end
    end
end %build the distance matrix for the neighborhood within 
% tmp = 0.05*(randn(size(ldist))); %added 5% error of normal distribution
% ldist=ldist.*(1+tmp);
% commented out the ablve two lines on May 8th


function [node_index]=find_neighbors(D,r,i,k)
%localDist(D,r) takes a distance matrix D and range r to generate for  
%node 'i' 0<i<size(D) a list that includes the nodes j such that D(i,j)<=kr. 
 
node_count=0;
for j=1:size(D)
    if D(i,j)<=k*r %look for all the nodes that is within distance k*r 
        node_count=node_count+1; %count the number of selected nodes
        node_index(node_count)=j; %save the original index of the selected node
    end
    %The saved node_index is in ascending order
end